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Abstract. We present an exact mapping between two simple spin models: the 
Fredrickson- Andersen (FA) model and a model of annihilating random walks with 
spontaneous creation from the vacuum, A + A 0. We discuss the geometric structure 
of the mapping and its consequences for symmetries of the models. Hence we are able 
to show that the upper critical dimension of the FA model is two, and that critical 
exponents are known exactly in all dimensions. These conclusions also generalise to a 
mapping between A + A <-> and the reaction-diffusion system in which the reactions 
are branching and coagulation, A + A ^ A. We discuss the relation of our analysis 
to earlier work, and explain why the models considered do not fall into the directed 
percolation universality class. 
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The one spin facilitated Fredrickson- Andersen (FA) model 5112] has been of interest 
recently as a simple model that exhibits dynamical heterogeneity IH ISl El IZ| • This 
model has a dynamical critical point at zero temperature that is not characterised by 
the divergence of any static lengthscale. As such, it represents a possible model of 
structural glasses, in which dynamical lengthscales seem to be large, but no large static 
lengthscales have been found 

An important advance was made by Whitelam et al. 4j. They showed that 
the master equation of the FA model can be cast in a field-theoretic formalism that 
resembles a simple reaction-diffusion system with branching and coagulation proceses. 
The dynamical critical point occurs when the density of the diffusing defects vanishes. 
The properties of the system near the fixed point can then be studied by the powerful 
methods of the renormalisation group. The analysis of Ref. jlj indicated that fluctuation 
effects are important below four dimensions and that the model is controlled by the 
directed percolation (DP) fixed point 0| between two and four dimensions. 

In this article we follow Ref. jl] in writing the FA model in a field theoretic 
language. The overall picture of a zero temperature dynamical fixed point that controls 
the low temperature scaling remains robust. However, we show that a somewhat 
hidden symmetry of the FA model means that the fixed point governing the scaling 
is Gaussian above two dimensions, and identical to that of annihilation-diffusion below 
two dimensions. Hence the DP fixed point is not relevant to the FA model. 

We determine the critical properties of the FA model by means of an exact mapping 
to a model of diffusing defects that appear in pairs from the vacuum, and annihilate in 
pairs: we refer to this as the AA (appear and annihilate) model. The mapping holds in 
all dimensions and at all temperatures. In one dimension the AA model is a (classical) 
Ising chain with particular single spin dynamics; above one dimension it is more familiar 
as the reaction-diffusion system A + A 0. The critical properties of the AA model 
were derived by Cardy and Tauber jH]; by using our mapping we can then apply this 
derivation to the FA model. The mapping also allows us to identify an (exact) duality 
symmetry of the FA model, which it inherits from the parity symmetry of the AA model. 
These symmetries are most simply expressed in terms of the master operators for the 
stochastic processes, which can also be interpreted as Hamiltonians for quantum spin 
models [El HH [12] . 

We will show that the mapping from FA to AA models is a specific case of a more 
general relationship between the reaction-diffusion processes in which the reactions are 
A + A ^ on the one hand and A + A ^ A on the other. (Note the presence of both 
forward and reverse processes in these reactions.) Our conclusion will be that these two 
systems are controlled by the same fixed point, and that their critical exponents are 
therefore identical. 

The form of the paper is as follows. We define the FA and AA models in Section ^ 
Both are models of hard core particles, but have generalisations in which particles can 
share the same lattice site; we also define these models (which we refer to as 'bosonic'). 
The mapping between the FA and AA models is determined by their symmetries: it 
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has quite a rich geometrical stucture, which we discuss in Section |21 In Section it is 
shown that the same mapping also connects the more general reaction-diffusion processes 
A+y4 ^ A and A+A ^ 0. An analagous mapping also exists between the corresponding 
bosonic models (Sectional). In Section HI we study the critical properties of the various 
models. We also verify the scaling in three dimensions, where the differences between 
our predictions and those of Ref. |3] are clearest, by numerical simulation. Finally, we 
summarise our results in Section El 



1. The models 

1.1. Model definitions 

In this section we introduce several models that we will consider in the remainder of the 
paper. We define them in terms of microscopic dynamical rules, before outlining the 
methods by which we represent their master equations and correlation functions. 

We begin with the one-spin facilitated Frederickson- Andersen model pj. This is 
referred to simply as 'the' FA model in what follows; more general FA models with 
facilitation by several spins, which exhibit more cooperative behaviour 0, will not be 
covered here. We define the model in terms of N binary variables, Ui G {0, 1}, on a 
hypercubic lattice in d dimensions. The Hamiltonian for the system is trivial: 

E = J2n^. (1) 

i 

We refer to a site with n = 1 either as an up-spin or as a 'defect'; sites with n = are 
thought of as down spins or 'empty'. The spins can flip with Metropolis rates if and 
only if at least one of their neighbours is in the up state. That is, for nearest neighbours 
i and j, we have 

Ojlj >■ Ijlj; rate c, , , 

Ijlj Ojlj, rate 1. 

Since we have detailed balance with respect to the Hamiltonian E, the (dimensionless) 
rate c parametrises the temperature according to c = e~^, with P = 1/T as usual. It 
also sets the equilibrium density 

Q 

{rii) FA,eq = rieq = -—■ (3) 
-L ~|~ C 

We have introduced the notation (■)FA,eq for an equilibrium dynamical average. 

To avoid confusion, we note that there are other versions of the FA model in the 
literature. In the original model jl] as described here the rate for flipping a spin i is 
proportional to the number of its up-spin neighbours. An alternative definition chooses 
the rate for flipping spin i to be independent of the number of neighbours in the up state, 
as long as there is at least one such neighbour jTHI- We do not believe that this choice 
makes any difference to the critical behaviour of the system, but our exact mappings 
apply only to the original model as defined above. 



Mappings between reaction- diffusion and kinetically constrained systems 



4 



We can also define a variant of tlie FA model in which the occupations Ui are not 
restricted to binary values, but may be any non-negative integers. We refer to this as 
the bosonic FA model, since the natural field theory for describing it has bosonic fields. 
This model has transition rates 

UiUj {ui + 1)%, rate cn^, , . 

{ui + l)nj UiUj, rate (n^ + l)nj. 

Here and throughout, parameters for the bosonic models are distinguished by tildes 
from those for the hard core ones. The bosonic model (jU obeys detailed balance with 
respect to the stationary state P({nj}) = Y\^e~^ c^"- /uil, so the stationary state density 
is 

(^i)FA,eq = C- (5) 

The bosonic stationary state is again a Gibbs distribution with Hamiltonian E and 
temperature defined by c = . However, in contrast to the hard core case it includes 
an a priori phase space weight factor l/nj, as appropriate for boson statistics. 

We will establish a mapping between the FA model and another model of diffusing 
defects. In this model, defects appear in pairs out of the empty state, and annihilate 
in pairs into it; they also diffuse freely across the lattice. We refer to this model as 
the AA model, since the defects appear and annihilate. Both the AA and FA models 
can be interpreted as reaction-diffusion processes. The AA model has explicit diffusion, 
combined with reversible annihilation A+A ^ 0. (We follow the standard nomenclature 
of reaction-diffusion models here, with A referring to our single species of particles, i.e. 
defects.) The model is defined for binary variables: for nearest neighbours i and j, 

IjOj Ojlj, rate 7c', 

Ijlj OiOj, rate 7, (6) 
OjOj — > Ijlj, rate 70'^, 

where we choose 



This model also obeys detailed balance with respect to the trivial Hamiltonian, E, at 
a temperature parametrised by c'. Note that, once a trivial overall scale for the rates 
has been removed, there are in principle two dimensionless rates, and the requirement 
of detailed balance with respect to E only fixes one of these, namely the ratio of the 
rates for appearance and annihilation. For now, we do not let the diffusion rate vary 
independently and instead tie 7 to c'. The mapping between the FA and AA models will 
then connect models with the same value of the single parameter c. The two-parameter 
generalisation of the AA model with an arbitrary diffusion constant also has a mapping 
to a generalised FA model, as we discuss in Section |221 For notational convenience we 
study the above 'standard' AA and FA models first in what follows. 
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Finally, there is also a bosonic variant of the AA model with rates 
UiUj — > (rij — + 1), rate 7c' n,, 

UiUj — > [ui — l)(r;,j — 1), rate jUiUj, (8) 
UiUj {ui + l){nj + 1), rate 70'^, 

where we take c' = c/A and 7 = 2. The stationary state densities in the AA models are 

c' 

Note that in the limit of small c, the relation ((7j) between c' and c becomes c' = c/4, 
and thus identical to the one for the bosonic models. More generally, the hard core and 
bosonic models should become effectively equivalent at low densities where multiple 
occupany of sites is very unlikely. We will exploit this correspondence frequently. 



1.2. Operator forms for the master operators 



It is convenient to write stochastic averages for systems such as those defined above 
in an operator formalism ^1]. This is a standard technique, so we largely restrict this 
section to definitions of the quantities that we will use later. 

When considering the bosonic versions of the FA and AA models we use a bosonic 



algebra with creation and annihilation operators on each site: [oj, aj] = [c 



0; 



6ij. 



ila]] 

The state {ui} is then associated with the vector Hil'^I)"' where 
|0) is the vacuum state which has all sites empty; the set of all 2^ state probabilities 
P{{ni},t) is mapped to the vector \ipit)) = J2{ni} ^{{''^i}^^) ni('^l)"'|0)- The individual 
probabilities can be retrieved via P{{ni}, t) = (0| Hil'^r/'^iOlV'l^))) and since they must 
sum to unity one has (e|^/'(t)) = 1 where 



(01 n- 



(10) 



is a 'projection state' that implements the sum over all possible system configurations. 
The master equation can then be written in operator form as dt\ijj{t)) = —C\ilj{t)), 
where C is known as the Liouvillian or simply the master operator. The off-diagonal 
elements — (0| nj('^rV''^iO'^ Iljl'^i)"' |0) S^^^ rates for transitions from state 

{ui} to {n'j^}, while the diagonal elements follow from the requirement {e\C = 0. Since 
the master equation is linear, it can be solved formally as \ip{t)) = e~'^*|'?/'(0)). If we 
specify the initial state as {rij} we can read off from this the probability of making a 
transition to state {n'j} in some time interval t: 



n 



n' 



-ct 



|0). 



(11) 



Expectation values over the stochastic dynamics can also be expressed in a simple form; 
for example, the average of some function f{{ni}) at time t becomes 

(/(KD) = E/(K})^(K},^) = (e1/(K})e-^*|^(0)), (12) 

{rii} 
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where fii = a\ai is particle number operator for site i. Similar expressions can be written 
for correlations functions involving two or more times, as illustated below. 

The hard core models with their binary occupation variables Ui have similar 



m a spin- 



relations but here the states are generated by operators sf and s^ e 
I algebra, with {sfY = = and sfs~ + s~sf = 1. The state vector is now 

m)) = Ein.}Pi{n^},t)U^istr\0), and conversely P{{n,},t) = (0| n.(C)"' I^W)- 
Conservation of probability requires (e|?/'(t)) = 1 with the appropriate projection state 
now being 

(e| = (0|n(l + ^n- (13) 



Transition probabilities and single-time averages take the forms 



n( 



n( 



|0), 



and 



(/(K})) = (e|/({n,})e-^*|^(0)). 



(14) 



(15) 



respectively, where now fii = sfs^. Occasionally it will be useful to write states and 
operators in notation analogous to Pauli matrices and spin vectors. Choosing a basis at 
each site as 



li). = |0). 
one has for example 




IT). = |i). 





Hi 




(16) 



(17) 



In principle one should write in this expression a direct product <S>jyti with identity 
operators at all other sites, but for ease of presentation we drop this here and below, 
along with site subscripts i where these are clear from the context. Our ordering of the 
basis states, while the reverse of the usual convention for spins, facilitates comparisons 
with other work on reaction-diffusion systems ^21- It also emphasises the analogy to the 
bosonic case, where the only natural ordering of the basis states is in order of increasing 
occupancy. 

It remains to give the forms of the master operator C for our models. Their matrix 
elements are easily derived from the relevant transition rates as explained above. One 
finds: 



^FA = Y\i4- ^)a]aj{ai - c) + [i ^ j] 



Uj, = (7/2) J2 M - i)^r(4 + 1)^74(^7 + c')^.n^7 - c) + {^^J 



(18) 
(19) 
(20) 
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£aa = (7/2) J2 [(4 - + + - c') + j)] , (21) 

{ij) 

where the sums run over all nearest neighbour pairs. The operators £fa and Caa for the 
hard core models have been written in a suggestive form that emphasises the connection 
with their bosonic counterparts £fa and £aa- 

Before leaving this section, wc note that the stationary states of our models have 
simple closed forms in the quantum formalism, viz. 

i i 

for the hard core and bosonic case, respectively. The latter is distinguished by a tilde 
as usual. Correlations in the stationary state then also take a rather simple form. For 
times ti, . . . ,tk that are in increasing order we have 

K(ti)ni,(t2) . . .ni,(tfe))FA,eq = (e|ni,e-^^^(*'=-*'=-i)ni,_^ • ■ ■ n,,e-^P^(*^-*i)n,Jc). (23) 

The A A model has an identical relation with £fa replaced by £aa and c replaced by 
c', and for the bosonic models one merely has to substitute for the master operator, 
projection and stationary state vectors their bosonic equivalents. 



2. Symmetries and mappings for hard core particles 
2.1. Detailed balance, parity and duality symmetries 

Having set up the operator formalism for dynamics, we now investigate some properties 
of the Master operators for these models. We first consider the effects of detailed balance, 
which tells us that the operator Ce~^^ is Hermitian (or more specifically symmetric, 
since all matrix elements are real). Here E = ^^fii is the (Hermitian) operator for the 
energy. Multiplying by ^^1'^ from the left and right shows that also 

H = e^^/2£e-^^/2 (24) 

is Hermitian. This is more useful than Ce~^^ since it is related to the LiouviUian C by 
a similarity transformation and so has the same eigenvalues. For the FA model we can 
write explicitly e~^^/^ = hi{c), where /ij(-) is the single site operator 



hAx) 



^^''^4^1 + «i «i 



1 

xV2 



(25) 



For the AA model we only need to replace c by d . The Hermitian forms of the 
Liouvillians are then 



fa 



'FA 



H 



AA 



c 



AA 



Their explicit forms make it evident that they are indeed Hermitian: for example 

^^FA = M - V^K4^j4{s7 -V~c) + {i^ j)] , 

(ij) 



(26) 



(27) 
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and -ffAA is similarly obtained from £aa in (|2Up by replacing the coefficients ±1 and ±c' 
by ±^/d. 

The above similarity transformation to a Hermitian form of the Liouvillians is 
convenient since it makes manifest the symmetries and conserved quantites of the 
systems. The mapping between FA and AA models relies on the fact that the 
Hermitian operators ifpA and Haa are related by the exact unitary (or, more specifically, 
orthogonal) transformation 

= U-'HaaU, (28) 



where 



1 VC 



VT+d L M y/T+d \ -VC 1 

with s^ = (s^ — s^)/22 as usual. Equation (|28p is the key relation from which most other 
results for the hard core models are derived; it is easy to verify by direct calculation. 
The operator U has a simple geometrical interpretation: it is just a rotation about the 
y-axis of the spin sphere, as illustrated in Figure Q below. 

From (j^Hj) we have directly a similarity transform between the corresponding master 
operators for the FA and AA models: 

CpA = V-'CaaV, (29) 




with 

= n^»' = ^'^2 ^ \l 'c' 

i 

We have exploited the freedom to introduce an arbitrary prefactor into Vi to ensure that 
both its columns add up to unity, i.e. 

((0|, + (l|,)i;, = (0|, + (1|,. (30) 

For the whole transformation V this implies that the projection state ()13|) is invariant 
under multiplication by either V or from the right, {e\V = {e\V~^ = (e|. So ()29p 
automatically maps a probability-preserving Liouvillian onto another one. 

Various relations between correlation functions in the two models can now be 
established. In addition to (jSOI) one uses the analogous property for application of 
Vi to the steady state vector on the right: 

|0). + c|l). ^ |0). + c^|l)^ , . 

1 + c 1 + c' ' ^ ^ 

and hence V\c) = \c') for the steady state ()22|) . For the simplest connected correlation 
function one has then, using the definition (jSHI) and the mapping (jSH)), 

= {e\{ni - neq)e"^P^*(nj - ncq)|c), 

= {e\V-^[Vini - r2eq)\/-^]e-^A^*[F(n, - n,^)V-^]V\c) . (32) 
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Equation (jSOI) implies that the leftmost factor of V~'^ can be absorbed into the projection 
state, while the rightmost factor V just changes c to c' in the steady state vector. Given 
that the number operators transform as 



-1 



VhiV 



one verifies also that 



-1 




((0|i + {l\i)vi{hi - n, 



:((0|, + (l|i)(ri,-n: 



eq/ 



and 



,i|0). + c'|l)i 



1 + C' 



|0). + c'|l)i 



so that overall 



-<q]K(0) -<q]> 



1 



(33) 



(34) 



(35) 



(36) 



X _)_ C ^'-'^'^"^ -eqjL-v\-/ -eqJ/AA.eq' 

It is then a trivial extension to show that for arbitrary connected stationary state 
correlation functions of a single time difference, 
I 



n, 



eqj 



r=l 



1 +C 



(/+-m)/2 




n, 



oq/ 



r=l 



(37) 



AA,eq 



However, a direct generalisation to stationary state correlations involving more than one 
time difference, or out-of-equilibrium quantities depending on more than one time, is 
not possible. This is because the transformation (jH^ of the number operator produces 
a non-diagonal operator which does not directly correspond to a physical observable. 
Only where the transformed operator is applied either to the projection state on the 
left, as in (jSH), or the steady state on the right, as in (j35|) . can such a link be made; 
otherwise more complicated relations result |T5]. 

The most useful aspect of the mapping ()29p is that it will enable us to reveal 
symmetries of the FA model which are 'inherited' from symmetries of the AA model. 
Specifically, it is clear from the dynamical rules of the AA model that the parity of the 
total number of particles in the system is conserved. Mathematically, we have that 



AA 



AA 



n2. 



(38) 



where -s^ = rij — | = °). Geometrically, the operator Hi 2sf simply produces a 
rotation of tt radians about the 2;-axis of the spin sphere. 

Since the FA and AA models are related by a similarity transformation, there must 
be a symmetry of the FA model that is equivalent to the AA parity symmetry. Applying 
the transformation to equation (jHHj) . we arrive at 

W-^CyaW, (39) 



FA 




Figure 1. The geometrical structure of the mapping U and the symmetries of the FA 

and A A models, in terms of the Hcrmitian operators Haa and i?FA- Points on the 
spin sphere represent states e"*"^/^ cos(0/2)|t) + e*'^/^ sin(^^/2)|4) where 6 and cf) are the 
usual polar and azimuthal angles. The black dots mark the position on the spin sphere 
of the zero eigenstates of the operators; these are 0j(|i)i ± ■\/c'|t)i) for the AA model, 
and 0j(|i)i + \/c|T)i) and 0- for the FA model. Since those states factorise over 
sites i, the figure can be read not just schematically as representing the entire A^-spin 
system, but also literally as showing the spin spheres for a single site. The rotation 
U is about the y-axis of the spin sphere (which points into the paper): it maps Haa 
onto HpA- The rotation n;(2sf) of tt radians about the z-axis maps Haa onto itself. 
Applying the mapping U gives the duality transformation W = U~^{Y[^ 2sf)U , which 
is a rotation of tt radians about the dashed axis and maps H-pA onto itself. In terms of 
the transformations in the main text, W is simply the image W after mapping Cya 
onto H-pA'- we show W here since its geometrical structure is simpler. 



with 

Note that = W, as expected for a symmetry deriving from the parity symmetry 

in the AA model. To understand more closely the effect of W note first that, in the 
AA case, the rotation Hj^sf maps the steady state vector \c') oc (3)j(|0)i + c'|l)j) to 
the vector (S)i(|0)i — c'|l)i) where the probabilities of all states {rii} containing an odd 
number of particles acquire a negative sign. The sum and difference of these two states 
gives the physical steady states for initial conditions containing an even and odd number 
of particles, respectively. In the FA case, W also maps two steady states onto each other: 
+ oc |0) = 0- \0)i. The symmetry thus links the 'conventional' steady 
state, which is reached for any nonzero initial number of particles, to the vacuum, i.e. 
the empty state; the latter is trivially a steady state since the kinetic constraints of 
the FA model forbid any transitions into or out of it. So while the original symmetry 
in the A A model connects steady states that are basically equivalent, with associated 
'domains of attraction' of equal size, the inherited symmetry of the FA model relates 
two very different steady states, with one having a domain of attraction containing all 
configurations except for the empty one. 
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The above relations between the FA and AA models, and their corresponding 
symmetries, can also be understood in terms of the associated Hermitian operators. 
They then have simple geometric interpretations, as shown in Figure ^ 

In the following we will continue to refer to 'the' steady state of the FA model 
as the one with nonzero particle density. The symmetry then allows us to relate 
the dynamics in this steady state to that in near empty configurations. This implies 
relations between the associated correlation functions. Proceeding as in (fH^. one has 
for example 



= {e\W~'^[W{ni - 



n. 



oq 



)w- 



[Win,-n,^)W-']W\c), 



{0\[W{h,-n,^)W 



[Win,-n,^)W~']\0). 



(40) 



Here we have used that application of W to the steady state vector |c) on the right gives 
a multiple of the vacuum state. The same is easily checked for the operation of on 
the projection state (e| on the left; the associated proportionality factors cancel because 
of overall normalisation. The transformation of the number operators is 



so that 



1 







{[Uiit) - ncq]K(0) - neq]) 



FA,eq 



(0| 



1 + c 



1 + c 



-1 

1-c 



-as, 



+ (1 



C Ui 



1 + c 



cs- 



1 + c 



|0). 



(41) 



Up to the overall numerical factor c/(l + c)^, the right hand side is of same form as p4|) : 
it is the probability of a transition between particular initial and final states, containing 
a single particle on sites j and i respectively. This relation generalises straightforwardly 
to correlation functions involving more than two spatial points: we have that 



n. 



oq; 



-1 



JJ(nj^(0) -rieq) 



FA,cq 



\l+m 



:i + c 



d+m, 



(o| 



n 



n 



Jr 



|0), 



(42) 



where the right hand side is again of the form (fT^ and gives the transition probability 
between an initial state with m particles and a final state with I particles. While this 
relation may not be familiar, it is closely related to the duality symmetry of the DP 
fixed point [S]. The latter is more usually expressed in terms of the dynamical action: 
see Section m 

In summary, we see that the transformation V maps the parity symmetry of the AA 
model onto an (exact) duality symmetry of the FA model. The mapping thus exposes 
a hidden symmetry which would not easily be recognised by looking at the FA model 
alone. 
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2.2. Models with additional diffusive processes 

The discussion so far has considered the FA and AA models, both defined in terms 
of a single parameter c. We now generalise our arguments to models with extra 
diffusive processes. This will show that our mapping applies more broadly between 
reaction-diffusion models with, respectively, reversible coagulation (i.e. coagulation and 
branching) and reversible annihilation (i.e. annihilation and appearance) processes. The 
generalised models will also allow us to elucidate the connection between our mapping 
and related earlier studies. 

Consider supplementing the FA model by an additional process 

Oilj hOj, rate D. (43) 

In particle language this is diffusion of a particle A to a vacant site, while the original 
processes are A ^ A + A (branching) and A+A^A (coagulation). This generalised 
FA model can therefore also be viewed as the reaction-diffusion model A + A ^ A. 
The new diffusion term in the master operator can be written as 

Cdm = D ^(s+ - sps-sfsptis; - sj), 

= -2DJ2 \s-s] + s\s] + sis] - 1/4] , (44) 

where s^ = (s+ + s~)/2 and, as before, = {sf — sf)/2i and si = hi — 1/2 = 
{sf s~ — s~ sf) /2. We recognise in (jUj) the Heisenberg model: see Ref. [TUj for a summary 
of the links between the properties of stochastic systems and their corresponding 
quantum spin Hamiltonians. 

For our purposes it is important to recognise that C^m has nonzero matrix elements 
only between states containing equal numbers of particles; it is therefore invariant under 
transformation with Hj^il^) ^ t\i&t the associated Hermitian operator is identical 
to £diff - and under the parity transformation J^- 2s^. Due to its Heisenberg form, £difT 
is also left invariant by any global spin rotation, and in particular by U . Combining 
these properties, invariance under V and W then also follow. Hence the structure of 
the preceding subsection is all preserved for the generalised models: the generalised FA 
model with diffusion rate branching rate c and coagulation rate 1 maps via V onto 
a generalised AA model with rates 

IjOj Ojlj, rate 7c' + D, 

Ijlj OjOj, rate 7, (45) 
OjOj — > Ijlj, rate 70'^, 

where 7 and c' depend only on c, as defined in ((Zj). We note that all generahsed FA 
models have conjugate AA models, but that AA models in which the rate for the diffusive 
process is less than 7c' cannot be mapped to FA models with positive rates. 

At this point, we make contact with two earlier studies. Krebs et al. jTT] studied 
the above generalised models at zero temperature but with nonzero D. Consistent with 
this, their mapping between the models is the limit of our mapping V for c, c' 0. 
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Henkel et al. ^21 implicitly had the full mapping V , but considered it only in the 
context of one-dimensional systems that are solvable by free fermions. The AA model 
then reduces to the Glauber-Ising chain (their model IV) and the relevant generalised 
FA model has D = 1 (their model II). Henkel et al. did not comment that the mapping 
applies to all dimensions and to arbitrary values of the diffusion constant. 

We illustrate the relation between our work and Refs. ^T] and in Figure |21 We 
parametrise the FA [A + A ^ A) models by the ratio of branching and coagulation 
rates c and the ratio of diffusion and coagulation rates D. In the AA case, appropriate 
dimensionless parameters are the ratio of appearance and annihilation rates c'^ and the 
ratio of diffusion and annhilation rates D' = c' + D/'j. All generalised FA models map 
onto generalised AA models with c' < 1 and D' > d . The standard FA model is D = 0, 
giving D' = d , while the pure coagulation/annihilation models of Ref. ^l] correspond to 
c = c' = 0. The free fermion condition of Ref. |I2j, finally, is the line D = 1 which maps 
onto the Glauber-Ising line D' = {1 + c'^)/2. (The Glauber-Ising chain has diffusion 
rate 1/2 and annihilation rate (1 + c'^)~^, giving the stated ratio.) 

To see the explicit link between our mapping and that of Ref. ^21? one notes that 
their free-fermion quantum Hamiltonian is directly in the form of the Liouvillian £aa for 
the AA model if its parameters are chosen as Di = D2 = 1, hi = h2 = {1 — c'^) /(I + c'^), 
?7i = 2/(1 +c'^) and 772 = riic'"^. Henkel et al. ^21 then show that the FA Liouvillian can 
be obtained by the similarity transformation Cfa = BCaaB~^, with B = Y[ibi and 

, _ / \ / cosh0 sinh0 \ f b \ 

'~ y 1/^ J y sinh0 cosh0 j yo 1/b j ' ^ ' 

Following through their analysis gives a = e^'^^'^c^^^'^, b = q-^'^/^c'^/^ and = i6, where 
tan 2^ = y/c. The latter condition can also be written as tan^ = v^, so that 9 is in 
fact the rotation angle associated with our mapping U. Inserting these values one has 
bi oc v^^ (with proportionality factor [(1 + c)/4]-'^/^) as expected by comparison with our 
Equation (j29j) . We note that Henkel et al. associate particles with down spins rather 
than up spins as here. However, because they also use the opposite (i.e., conventional) 
ordering of the two local basis vectors ||) and ||), the matrix representations of all 
operators are the same. 

To summarise, we showed in this section that the FA and AA models have quite 
a rich geometric structure underlying their symmetries and the relations between 
them. These relations further extend to a general mapping between reaction-diffusion 
models with coagulation and branching {A + A ^ A) and annihilation and appearance 
[A + A ^ 0). We expect the critical behaviour (at small particle densities, i.e. 
low temperature) of these models to be determined by their symmetry properties. 
However, the hard core constraint that allows only one particle per site makes an 
explicit renormalisation group analysis of such critical properties awkward. We therefore 
show next that the bosonic models, where this constraint is removed, have analogous 
symmetries and mappings between coagulation and annihilation models. 
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Figure 2. Sketch of the mapping between the two-parameter famihes of generahsed 
FA {A + A ^ A) and AA {A + A ^ 0) models. Generahsed FA models (with non- 
negative c) map onto generalised A A models with < c' < 1 and c' < D' . The 
standard FA models with D = map to models on the standard AA line c' = D' . 
The lines c — Q and c' — Q correspond to diffusion-limited pair annihilation (DLPA) 
and diffusion- limited pair coagulation (DLPC) respectively (llj . In one dimension the 
models arc solvable by free fermions on the lines marked 'ff'. These lines are given by 
D = \ and D' = {1 + c'^)/2, and the mapping transforms them into each other |12|. 
In the region marked MU the mapping is unphysical: such generalised AA models do 
not have FA counterparts with positive rates, though it seems unlikely that this would 
have physical consequences for the behaviour of the corresponding AA model. 



3. Bosonic models 



The bosonic models introduced in Section[T]have similar properties to those discussed for 
hard core (spin) models in the previous section. The Liouvillians again have Hermitian 
analogues defined hj H = e^^/^£e~^^/^. The energy operator is now E = ajflj, so if 
we define the bosonic version of hi(x) as 



then 



FA 



FA 



H 



AA 



c 



FA 



(47) 



(4^ 



HpA has the same form as (fT^ except that the numerical constants 1 and c are replaced 
by Vd, i.e. 



ai-Vc) + {i ^ j)]> 



(49) 



and -f^AA is obtained similarly from (j2H) . This is analogous to the hard core case, but 
easier to see for the bosonic models since the transformation by hi{x) simply rescales 
particle creation and annihilation operators, according to the first of the relations 



e-- ^F(a, a^)e-^'"' = F{e-^a, e^a\ , 
e^"^""^ F{a, at)e-^"-^"' = F{a - /i, + A). 



(50) 
(51) 
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One expects that at low particle densities Ueq = c/(l + c) ~ c the constraint of at 
most single occupancy in the hard core FA model will be irrelevant, so that it becomes 
equivalent to the corresponding bosonic model with c = c; the same argument applies 
to the hard core and bosonic AA models. This physical reasoning 4] can be further 
supported by a large-S* expansion of the hard core models ( |Appendix A ). 

We now discuss the mappings between the bosonic models and their symmetries. 
The main conclusion is that the structure of the hard core models is preserved in their 
bosonic counterparts. The basic mapping between the two bosonic Hermitian operators 
is 

^FA = U-'HaaU, (52) 
with the unitary operator 

i 

This is easy to verify, bearing in mind that c' = c/4 and 7 = 2: from (j3T|) . the 
transformation by U shifts all and a] in -^aa by —y/d/2 = —\f^. If one uses a 
basis of bosonic coherent states then the mapping is a translation in the complex plane 
that parametrises these states. This is the analogue of the rotation of the spin sphere 
generated by Wj, consistent with the intuition ( |Appendix A ) that the bosonic models 
effectively 'flatten' the spin sphere onto the complex plane of coherent states. 
Combining ui and hi we have a relation between the Liouvillians, 

£fa = y-^^KPy. y = \{v,, v, = kCc'Whr^-c) = 2-'^-^e^^-^4)/2^ (54) 

i 

where the explicit form of Vi follows using ()50|) and ()51|). These relations also show 
that the transformation by V is simply a combined shift and rescaling of the bosonic 
operators 

V-'a,V = i (a, - 1^ , V-'alV = 2 (al - . (55) 

The mapping between the two Liouvillians again relates the parity symmetry of the 
bosonic AA model, 

i:AA = i-l)^^''''''jCAAi-l)^^''^^\ (56) 

to the duality symmetry of the bosonic FA model 

Cfa = W~^CfaW, = t/-i(-l)S.4«>V^. (57) 

One finds explicitly 

W = W-^ = Y[{-iy^^^e^^-^''l (58) 

i 

and the duality symmetry transforms the bosonic operators as 

W-^aiW = c-ai, W-^alW = 1 - a|. (59) 
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The mappings V and W allow us to establish the analogues of (jHTj) and (ji^ . The 
first of these relates steady state correlations in the bosonic FA and AA models via 



FA,eq 



i-r=l 



c') 



l-r=l 



(60) 



AA,eq 



Here the v must label sites that are all distinct from each other, as do the jr, though 
the two sets may contain sites in common with each other. The prefactor on the right 
agrees with the low density hmit of the one in ^7\i . supporting our intuition about the 
equivalence of hard core and bosonic models in this regime. 
The duality symmetry of the bosonic FA model results in 

- / -I r- m 



r=l 



FA,eq 



-1 



■{0| 



n 



-C-p\t 



n 



Jr 



|0). 



(61) 



This again relates steady state correlations to transition probabilities between specific 
initial and final states; the prefactor approaches the one in (jl^ for c = c — > 0. 

Following our discussion in Section lT^ one expects that the structure of the bosonic 
mapping will be preserved if an extra diffusive process is added to both models. This is 
easily verified. It is also immediate to show that the mapping is unchanged if we add 
on-site branching and coagulation processes, as long as we retain detailed balance for 
the whole model. These processes then map to on-site appearance and annihilation in 
the generalised AA model. (Recall that our standard bosonic FA and AA models have 
processes that always act on pairs of sites.) 

The similarities between Eqs (147116111 and ()25IIH9j) show clearly that the bosonic 
models have the same structure as those with hard core exclusion. In the next section 
we consider the critical properties of these bosonic models; from our arguments the 
critical properties of the hard core models should be identical, and we check this by 
comparing our predictions to numerical simulation. 



4. Critical properties 

We have shown that the FA and AA models are linked by an exact mapping. Now, both 
models have scaling behaviour at small defect densities that is characterised by the fixed 
point of a renormalisation group (RG) fiow. In Section 14.11 we use the mapping of the 
previous section together with known results to show that the FA model has upper 
critical dimension dc = 2; this conclusion also applies to the generalised FA model, i.e. 
the react ion- diffusion model A + A ^ A. The critical scaling is then characterised by the 
well-known mean-field (Gaussian) exponents in d > 2; we also derive exact exponents 
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below dc that coincide with known results in one dimension. Our results differ from 
earlier studies in two and three dimensions: in Section we therefore use simulations 
to confirm the predicted mean-field scaling in d = 3. In Section 14 .HI we derive some 
analytical results for the scaling limit of correlation functions in > 2. Finally we 
discuss in Section 14.41 the scaling of the persistence function since data for this were 
used in Ref. j3] to support the argument that non-mean-field fluctuation corrections are 
significant in three dimensions. 



Renormalisation group analysis 

The critical properties of the bosonic AA model were established by Cardy and 
Tauber [Oj: in their notation the model corresponds to /c = 2, r > 0, 0"™ = 0. We 
write the generating functional for dynamical correlations in the stationary state as a 
path integral on the lattice 



^AA 



where (fu and (p\^ are time-dependent conjugate fields at each site i. The 'Lagrangian' 
LaaH'^u}, {vlt}] is obtained from Caa by replacing — > (pu and aj (fl^. and depends 
on all the fields at a single time. 

Taking the continuum limit, the lattice fields {(pu} are promoted to a field (pxt 
depending on spatial position x and time t. The generating functional becomes 

I?[0,t,0ye-^^^[^-'^-l (63) 



^AA 



where the functional ^aa is known as the dynamical action. Including gradient terms 
up to second order gives 



+Ao(0l, - - p')(l + ^oVV2)(0L + + P) (64) 

where we have neglected boundary terms; p' is the steady state density (proportional to 
c'), /o is the microscopic lengthscale (lattice spacing) and Aq is a bare coupling constant 
that sets the microscopic timescale. The dimensions of Aq are [time]~^ [length]'^; the field 
0^ is chosen to be dimensionless and (j) has dimension of [length] 

While the above factorised form for the action was useful for the exact mappings of 
the previous sections, the RG calculation requires us to separate the terms in the action 
that correspond to different physical processes. We write 

Saa[<PA^] = j ^'xdt |0t^(9i-Aop7^v2)0.t + Ao[(0L)'-l](0'.-p'') + ^AA,i} (65) 
where 

^AA,i = (Ao/^/2) [4,0,, V^^L^x. + ( V0,,)2 + {p'Vct>lf] (66) 

(we continue to ignore boundary terms when integrating by parts over spatial degrees 
of freedom). Physically, we recognise the first term in (|65|) as a diffusive propagator 
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for the excitations and the second term as local appearance and annihilation processes. 
The terms contained in Laa,i will be irrelevant in the RG sense since their only effect is 
to enforce the fact that appearance and annihilation of excitations take place on pairs 
of adjacent sites and not on single sites. (In terms of the RG calculation these terms 
modify the spatial structure of terms in the action that are already present, but they 
are not responsible for new terms, or for any singular behaviour.) We therefore neglect 
Laa,i and arrive at the action considered in Ref. P for the case {k = 2, r > 0, am = 0) 
described above. 

We now follow Ref. in renormalising this dynamical action. The parity symmetry 
of the bosonic AA model is 

SAA[(pxt, (pit] = SAA[-<i)xt, -<plt], (67) 

and must not be obscured by making any shift of the fields ;9]. The symmetry is clearly 
preserved under the RG flow so only terms in the action with this symmetry need be 
considered. Power counting then shows that the upper critical dimension will be two. 
Above d = 2, therefore, the critical exponents have their mean-field values 

(z,i/,/5),>2 = (2,1/2,1). (68) 

Here we have defined P by the scaling of the steady-state density \imt^oo{ni(t)) ~ p'^, 
and z/ by the correlation length scaling ^ ~ p'"- (The notation a h means that a is 
proportional to h in the scaling limit, i.e. close to the critical point at p' = 0.) Note 
that in Ref. P these exponents were defined in terms of the control parameter r ~ p'^ 
and thus differ by factors of two from ours. Our convention is more appropriate for 
comparison with the FA model, where the steady-state particle density is the natural 
control parameter. We also note that the free propagator and hence the bare diffusion 
constant Dq (the constant multiplying the term (j)^'V'^(f) the action) depend explicitly 
on p'. In the usual RG analysis Dq is set to unity, so we define the exponent z via the 
scaling of typical relaxation timescales r measured in units of -Dq ^, 

Dqt ~ ~ p'-'". (69) 

The scaling of the times r in absolute units is then governed by an exponent different 
from z: Dq (x p' gives 

r ~ ^7^0 ~ p'-^^'r (70) 

We next show that below d = 2 the exponents are exactly 

(z,i/,/5),<2 = (2,lMl). (71) 

This is consistent with the exact scaling in one dimension and also with the naive scaling 
estimate r ~ ^-1-2/^ fQj, ^y^q -p\ model 0. In Ref. 0, u and P were given only to first 
order in a loop expansion: in our notation the results were 

z = 2, p = du, u = 2/yr = l/d + C(2 - df, (72) 

for d < 2. Here y^. = 2d + 0{2 — dY is the scaling dimension of the coupling r ~ p'^. 
(The relation (3 = dv comes from the scaling relation (5 = zv/a where a is the exponent 
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for the decay of the density after a quench to criticahty, which equals d/2 fU].) However, 
we can determine exactly since detailed balance fixes the steady-state density to the 
value c' oc p', so that (3 = 1. Thus we have the exact result z/ = 1/d, or 

Vr = 2d, (73) 

and the exponents (ffTjl are exact. Equation ()73j) can be confirmed by summing the 
geometric series of loop corrections that contribute to the scaling dimension y^-. 
Proceeding to the bosonic FA model, the (unshifted) dynamical action is |3] 

Sfa = j d''xdtct>lM.t + /Wo(0L - - P)(l + /^VV2)0L0xi, (74) 

where p oc c is the steady state density and po is the bare coupling (with dimension 
[time] "^[length]''); Iq is the microscopic lengthscale as before. The property of detailed 
balance manifests itself as an invariance of the action, 

5'FA[0xt,0L] = S'FA[p0l_t,P"Vx-t]- (75) 

This of course has an analogue in our earlier operator notation, where we recognised 
detailed balance as the fact that £FAexp(— /3^) is Hermitian. Since exp(— /5^) = 
gEi«iii^ this implies exp(— ^fa*!^) = c"^*"^'"' exp(— £p^t)c^»'^»"* which is the promised 
analogue of (f75j) [recall (jM^j) and note that /"'"(a, a"'') = /(a"l',a)]. 

The duality symmetry as it is normally stated for systems with a DP fixed point 
(but not necessarily with detailed balance) is |H] 

Sya[<Pxu (pit] = Sfa[p{1 - 0x,-t)> 1 - P"Vx,-i]. (76) 
Like detailed balance this involves time reversal; in terms of the Liouvillian it relates £fa 
to its conjugate £p^. To arrive at the duality mapping for the FA model, we combine 
the preceding two symmetries of the action to obtain 

SFA[(pxt, (pit] = Sfa[p - (px,u 1 - (pi,t\- (77) 

In terms of operators, this symmetry now relates £fa directly to itself, without any 
conjugation. A comparison with reveals that it corresponds directly to our earlier 
transformation W from (jFTj) . 

Both detailed balance and duality symmetries are preserved under renormalisation 
if all terms in the action are retained. It is then crucial to follow Ref. [0] in choosing 
a basis for the RG equations that reflects this fact. The solution is to make the 
transformation defined by V in the previous section and to write as in (j55|) 



Hence we can write the dynamical action in terms of these new fields: 
^fa[^,V^^] = J d^xdt^lAi^^t 

+fio{^l - - P/4)(1 + llV'mi^l + 1)^,, + p/4) (79) 

This is identical in form to Saa as given in Equation (jMj) . and so the FA model 
renormalises exactly as the AA model, consistent with their correlation functions 
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obeying the simple relation (jUUI) . The significance of the basis used is that the resulting 
RG equations respect both the detailed balance and duality symmetries of the FA 
model. If one makes any transform other than (|7H|l . using for example the standard 
shift = 0"'' — 1 as in Ref. then these symmetries are obscured and one is led to the 
conclusion that the FA model is controlled by the DP fixed point between two and four 
dimensions. 

4-2. Simulations showing Gaussian scaling in d = 3 

We now confirm our above exact predictions of mean-field (Gaussian) critical exponents 
in d > dc = 2 hj equilibrium simulations of the FA model in d = 3. We explain 
how the duality relation PT|) together with an appropiately adapted continuous time 
Monte Carlo (MC) algorithm allow us to probe citical properties well beyond the regime 
accessed in previous work. A comparison of our data with analytical scaling forms is 
also given. 

A convenient observable for extracting the relaxation time scaling is the normalised 
two-point susceptibility X2{t) = {AE{t)AE{0))/{AE{Oy) where AE = E-{E) denotes 
the fiuctuation of the energy away from its equilbrium value. Since we only consider 
the stationary state dynamics of the FA model in this section, we drop the subscript 
'FA,eq' on the averages. Substituting the simple form (P) of the energy function and 
using translational invariance, X2it) can be recast in the form 



with j an arbitrary reference site and rieq = (ni), Equation (jH)). The second equality in 
(jHUj) follows from the duality relation pij) . It shows that the stationary state average 
defining X2{t) has a dual counterpart in the dynamics near the empty state. This is a 
tremendously useful fact: instead of having to simulate an equilibrium system containing 
hundreds or even thousands of defects we simply initialise with a single defect at site j. 
According to (jHOI), X2(t) is then given by the probability that this state evolves under 
FA dynamics into one containing again a single defect (at an arbitrary site i), that is 
any configuration with E = 1. 

To measure with similar efficiency a dynamically growing lengthscale in the FA 
model one can consider the mean squared displacement associated with two-point 
correlations, 



Here x{i) denotes the position vector of site i, and j is again an arbitrary reference site; 
we set the lattice constant to unity so that x{i) e Z"^. Note that due to normalisation 
this lengthscale can exceed the equilibrium dynamical correlation length, which would 




(80) 




(81) 



conventionally be extracted from the maximum of X2{t)r'^{t). As in ()80|1 we apply 
the duality PT|1 to map the stationary state average in (jHTjl onto {0\sl'e~'~^'^^s^\0) , 
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and perform simulations of the dual problem rather than the time-consuming direct 
equilibrium simulations. 

The duality mapping increases the efficiency of an MC algorithm to such an 
extent that the limiting constraint in practical simulations is system size rather than 
computational speed. Although we initialise with just a single defect, its trajectory 
under FA dynamics explores the simulation box and must not be biased by finite size 
effects. For conventional, lattice based algorithms the required system size N = L'^ 
can then quickly exhaust the available memory of standard computers (say 1Gb). To 
overcome this problem we used a coordinate-based variant of MC. Instead of storing the 
occupation numbers Ui of lattice sites i we keep track of the actual coordinates Xa 
of each defect a in the 'virtual' simulation box {1, 2, ... , L}'^. The memory efficiency of 
this approach against a lattice based code is 0{M/N) if there are M defects. Hence it 
is useful for simulations of problems with low defect concentration, for instance, direct 
equilibrium simulations at low Ugq or the dual dynamics near the empty state. In analogy 
to lattice based codes [inillZI one can set up a 'reverse lookup scheme' where continuous 
time MC steps have 0{1) computational complexity. This is accomplished by storing 
the list of defect coordinates {xa}aLi in a hash table. So the problem of whether any site 
x{i) in the virtual simulation box is occupied by a defect can be decided in 0{1) time, 
just as for a lattice based code. This coordinate-based continuous-time MC approach 
is memory efficient while yielding computational speeds comparable to a lattice based 
code. This allows us to exploit fully the simplifications afforded by (jHUj) and (jHTj) . 

Simulation results for branching rates c = 10~^, 10~^, 10~^ and 10~^ are shown in 
Figures El and ^ below. We have used a formally infinite virtual simulation box so that 
via (j8n|) and (j8T| we are probing the equilibrium dynamics of the FA model in the 
thermodynamic limit. In other words, our results are guaranteed to be free of finite size 
effects and can be compared directly to the scaling predictions and (f7n|) . For this 
comparison we recall that at small defect densities the bosonic and hard core FA models 
have similar behaviour, so that p oc c ~ c is proportional to c. Since we know that 
(3=1 exactly, it is sufficient to show that z = 2 and zu = 1 to demonstrate Gaussian 
scaling. We first verify the dynamical critical exponent z = 2 which implies for the 
growing length scale r^(t), 

r^(t) ~ {Doty/^ ~ {ay/'. (82) 

For all values of c considered the data in Figure El show a linear dependence of r'^{t) on 
at, over up to eight decades in at; in fact our data are fully consistent with r^(t) = 2dDot 
where d = 3 and Dq = c/2. This value of Dq confirms the expected scaling Dq ~ c; 
we return below to a derivation of the prefactor in Dq = c/2. It should be emphasised 
that the kind of spatio-temporal scaling used here determines the dynamical exponent 
directly from its definition; measuring the ratio between exponents for the correlation 
length and relaxation time is a valid procedure only if the diffusion constant is slowly 
varying at criticality. 

Consider next the simulation data for X2{t) shown in the inset of Figure HI Scaling 




Figure 3. Simulation data (symbols) for equilibrium dynamics in the d = 3 FA 
model with branching rates c = 10"'^, 10~^, 10"'^ and 10~^ obtained from a coordinate- 
based continuous-time MC algorithm applied to the dual form of r'^it), Equation H81|l. 
Results are averaged over 10^, 10^, 10^ and 5x lO"* samples, respectively. The dynamical 
lengthscale increases according to a diffusive law r'^{t) = 2dDQt and with diffusion 
constant Dq — c/2 (full line and symbols). Error bars are significantly smaller than the 
symbol-size except where they are shown explicitly (data for c = 10^^ with ct < 10~^). 
The dashed line represents DP scaling of the dynamical lengthscale as discussed in the 
main text; it is inconsistent with the data. 



arguments predict that X2{t) should be a universal function of t/r, with r the relaxation 
time. From ()70|) we expect for Gaussian exponents 

r(c) ~ c-^-'" = c-\ (83) 

Collapse of data for X2{i) under this prescription is shown in Figure lU While there are 
still noticeable pre-asymptotic contributions at c = 10~^ the data for c = 10~^ and in 
particular c = 10~^, 10~^ seem to have converged to the final scaling form of X2{t) to 
within our numerical accuracy. To confirm the critical scaling in more detail we next 
give a theoretical analysis that predicts the precise shape of Xiif) in the critical limit 
0. 

4-3. Scaling analysis for two-point functions in d > 2 

We now show that the c — > limit of the equilibrium two-point susceptibility for FA 
models in ci > 2 dimensions is 

X2it) = exp{-Kdch), (84) 




Figure 4. Simulation data obtained from the dual representation (|8()|l of the two- 
point susceptibility X2 (t) and with the same parameters as in Figure 01 The inset 
shows X2{t) versus t. At each c data are collected up to times t = 0.5/c^, that is up 
to t = 0.5 X 10^^ for the lowest c. The main plot demonstrates data collapse with the 
Gaussian scaling of the relaxation time r ~ c^^. The full line represents our analytical 
prediction for the c ^ scaling of X2{t), Equation (|84|l : there are no fit parameters. 
Error bars are largest for the data at c = 10^®; even there we estimate a relative error 
below 1%. 



where = d [Pd{2, 0, . . . , 0) + 2(rf - l)Vd{l, 1, 0, . . . , 0)], and Pd(x) is the probabihty 
that a pair of random walkers with initial separation x will never meet: this function is 
considered in [Appendix B Equation ()84j) is consistent in particular with the scaling 
expectation that X2{t) should be a function of t/r with r ~ c^^. In d = 3 and 
substituting the survival probabilities ()B.6jl and ()B.7|1 one has K3 ~ 5.80961. The exact 
c ^ scaling form ()84j) for this dimensionality is shown in Figure 0] and is in perfect 
agreement with our simulation data. The remainder of this subsection is devoted to the 
derivation of ()84|) : we return to more general discussion in Section 

A systematic analysis of the dynamics in the critical limit c ^ requires a careful 
distinction between the timescales involved. We will need unsealed time t as well as 
the scaled time variables x = ct and y = c^t whose 0{1) increments correspond to time 
intervals 6t of 0{c~^) and (9(c~^), respectively. For time intervals 6t = 0{1) the limit 
c — > is clear: the rate c for branching processes 1^0^ — > Ijlfc, where is a nearest 
neighbour (NN) of j, then vanishes. Only the coagulation processes Ijl^ IjOk or 
Ijlfc Ojlfc can then take place, each occuring with rate unity, and the lifetime of 
excited states like 1^1^ that contain two (or more) NN defects is C(l). 
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The dynamics over intervals with 6x = 0{1), i.e. St = 6x/c = 0{c~^), is rather 
different. Processes involving branching events are then possible. However, after such 
an initial event the rate for reverting to the original state is 0{1), while the rate 
for an additional branching event is 0{c). This implies that the probability for the 
latter to occur first is 0{c). Since the rate for the first branching event is already 
0{c), the effective rate for a sequence of two branching events is then O(c^). On 
our 0{c~^) timescales this type of process can be neglected, and we need only concern 
ourselves with processes involving a single branching event. One possible process is then 
IjOk Ijlfc Ojlfc, with j, k NNs. Its rate is c x i, c being the branching rate and 
I the probability for the particular, subsequent coagulation (rather than Ijl^ IjOfc, 
which would take us back to the initial state). The probability for this process to occur 
during the interval St is therefore ^c6t = ^6x. Note that as the lifetime of excited states 
is 0{1) as argued above, we have vanishing probability 0{c) of finding the system in 
an excited state 1^1^ at any given moment in time. Therefore, in the limit c — 0, the 
intermediate excited step of our process becomes invisible and we have effective diffusion 
IjOk Ojlk with probability ^6x. 

The second type of process allowed on the 0{c^^) timescale is the excitation of a 
defect at a site that has more than one NN defect, followed by a cascade of coagulation 
events leading to another configuration where all defects are isolated. Consider first the 
case where the new defect, at site k (say), has two NN defects at sites j and /, so that our 
process leads from Ijlfcl^ to IjOfcO; or OjlfcO; or OjOkli- The corresponding probabilities 
are ^Sx, ^6x, ^6x, respectively: after the initial excitation (probability 2c6t = 26x) at 
site k, there is a probability of 1/2 that the defect at site k, which has twice the down- 
flip rate of those at j and /, will not flip down first. (If it does, we have returned to the 
original configuration and can ignore the process.) There is then probability 1/2 that 
j will flip before /, and probability 1/2 for each of the remaining defects to flip first, 
resulting in the overall probabilities given above. 

To summarise the dynamics on the 0{c~^) timescale, we have diffusion 1^0^ — > 0^1^ 
with rate | per interval of rescaled time x. The processes 1^-0^1; Ijlfel; — > . . . 
discussed above can then be represented consistently as produced by a diffusion event 
followed by (in the limit c — *• 0) instantaneous coagulation. For instance, IjOfcl; — *• IjOfcO; 
amounts to diffusion 0^1; 1^0; (with rate |) followed by coagulation 1^1^ IjOk 
(which has probability | of occurring before 1^1^ 0^1^), giving the overall rate 
i X I = i obtained above. One can check that the same decomposition into diffusion 
and instantaneous coagulation also holds for processes involving an initial excitation at 
a site with more than two NN defects. 

Our discussion so far leads to the following conclusions. Clearly any single-defect 
state s^|0) is blocked on the 0{1) timescale. Therefore X2{t) = 1 and r^(t) = in 
equilibrium and for C(l) times, according to the dual representations (jHn|l and (jHTjl . 
But on the 0{c~^) timescale the defect diffuses away from its starting site j. This 
still predicts X2(t) = 1 since diffusion conserves the number of defects E = 1. The 
quantity r^(t) in dual representation, on the other hand, now measures the mean squared 



Mappings between reaction- diffusion and kinetically constrained systems 



25 



displacement of a single, diffusing defect which is given by r'^if) = 2dDQt. Our result 
that the diffusion rate is ^ in time units of x = at tells us that Do = c/2. This is the 
exact low-density scaling of Dq and precisely what we found in our simulations; compare 
Figure 01 

Let us now turn to dynamics on the (9(c~^) timescale, corresponding to 6y = c^5t = 
(9(1). During such a time interval there is now a nonzero probability for the occurence 
of processes involving two successive branching events, where one defect is first excited 
next to an existing one and another defect is then created on a NN site to either of 
the other two. In order to determine the fate of this defect triple Ijlfclz we consider its 
evolution on the faster 0{c~^) timescale. This argument is analogous to the one above, 
where a defect pair is created on the 0{c~^) timescale but we have to look at 0{1) 
times to determine its relaxation. In further analogy we note that an arbitrarily small 
but nonzero increment 6y on the (9(c~^) timescale corresponds to an infinite increment 
6x = Sy/c of 0{c~^) time as c — 0. The following possibilities then arise: because NN 
defects coagulate instantaneously on the 0(c~^) timescale there is a probability of | for 
immediate relaxation Ijlfcl/ Im where m = j or k or I. Subsequently the defect 
diffuses with diffusion rate | for a time 5x — oo, moving arbitrarily far from its initial 
position. With the remaining probability of | an instantaneous relaxation of the middle 
defect takes place, Ijlfcl; — ^ Ij^kh- The defects Ij and 1; can now diffuse independently 
for an effectively infinite interval 5x of (9(c~^)-time; if they do not coagulate in the 
process, their distance grows without bound and we can say the original defect has 
branched irreversibly. We derive in Appendix B the 'survival probability' Vd{x) for this 
outcome; here x = x{l) — x{j) is the initial separation of the diffusing defects. 

We can now assemble the probability X5t that during a time interval 6t = 5y/c^ a 
single defect irreversibly branches into two defects. Starting from a single defect there is 
a rate c for branching on a given NN site. Since (on a hypercubic lattice in d dimensions) 
there are 2d such states, the overall branching rate is 2dc. The probability for a second 
branching event on a neighbouring site to take place before either of the two possible 
relaxations back to a single defect is c/2. A cluster of two NN defects has Ad— 2 NN sites, 
two of which lead to a hnear and 4((i — 1) to an angled defect triple Ijlfcl;. Altogether 
the rate for creation of a linear triple is 2d(? while it is 4:d{d — l)c^ for an angled one. 
In either case, we need to multiply by the probability | of the middle defect k relaxing 
first, leading to a pair of next nearest neighbour (NNN) defects. In terms of the eventual 
survival probabilities Vd{x) of this pair we thus obtain 

X = d [Vd{2, 0, . . . , 0) + 2(t/ - l)Vd{l, 1,0,..., 0)] c\ (85) 

This is a nontrivial result. To clarify its intuitive meaning, note that we are considering 
initially a single defect s^|0). The following trajectories are then possible during a time 
interval 5t = Sx/c = Sy/c^ on the 0{c~^) scale: (i) no branching occurs. The defect 
diffuses with rate | for a time 6x = Sy/c ^ oo in the c — > limit. The defect at 
the beginning and end of the interval 6t are then completely decorrelated. (ii) With 
probability d{2d — l)c^5t = d{2d — l)6y the defect branches into a pair of NNN defects. 
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Again, on the 0{c~^) timescale this pair has an infinite time 6x available to diffuse 
through the system, (ii.a) There is a finite probability that the pair coagulates during 
this diffusive motion. In this case we have within the time interval 6t a 'bubble' in the 
space-time diagram of the defect trajectories where the initial defect separates into 
two but these re-coagulate shortly afterwards. The temporal extent of this bubble, i.e. 
the time during which the two defects exist, is 0(c~^). The probability of detecting 
such a bubble on the (9(c~^) time scale is therefore vanishingly small in the limit c ^ 0. 
This is in analogy to excited states becoming invisible on the 0{c~^) timescale. Hence 
the trajectories (i) and (ii.a) cannot be distinguished on the 0{c~~^) timescale. (ii.b) 
The defects may diffuse forever {6x —>■ oo) without encountering each other; this means 
that the trajectory branches irreversibly on the 0{c~'^) time scale. Due to the existence 
of bubbles the probability for this event, A5t, is renormalised relative to the 'bare' 
probability d{2d — l)c^5t for an initial branching event. Each defect in the resulting 
pair has travelled an infinite distance during 5t and so completely decorrelates from the 
initial defect. 

We can now make predictions for the dynamics on the (9(c^^) timescale. First we 
note that in c? = 1, 2 the survival probabilities Vd{x) vanish, see [Appendix B Hence 
A = and defect trajectories do not branch on the 0{c~^) timescale. However, in c? > 2 
the Vd{x) are finite, and hence so is Hd = A/c^. Now X2(t) is just the probability that 
the number of defects has not increased during the time interval t, which means that no 
irreversible branching processes have taken place. The rate for occurrence of the latter 
being A, it follows that X2(t) = exp(— At) = exp(— K^c^t). This completes our derivation 
of Equation 

One can go further and extract the probability Pt{E) of having ii^ > 1 at time t. 
Since each irreversible branching event produces infinitely separated defects, these will 
then continue to branch independently in the same way. Thus, if E defects are present, 
the rate for generating an additional one by irreversible branching is EX. This gives the 
master equation 

dtPtiE) = iE-l) XptiE - 1) -EXptiE). (86) 

This can be solved straightforwardly, for example by ^-transform, to give 

p,iE) = - e-^'f-\ (87) 

The number of defects thus has an exponential distribution; the most likely outcome 
is = 1 at any time and has probability Pf(l) = X2{t) = e~^*. The average number 
of defects, on the other hand, grows exponentially as {E)t = e"*"^*. We have checked in 
our simulations the full form of Pt{E) (data not shown), and found excellent agreement. 
The exponential increase in the number of defects limits the time range that can be 
conveniently simulated; the lowest value of X2{t) that can be measured reliably is of the 
order of the inverse of the maximum number of defects one is prepared to track. 
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4.4- Persistence functions and discussion of earlier data 

In order to characterise a dynamical fixed point in the presence of a time-reversal 
symmetry, we must determine three independent exponents. One such set is {z,zh',P). 
Our simulations as well as the scaling analysis have confirmed the RG prediction 
that z = 2 and zu = 1; further, /3 = 1 is known rigorously from detailed balance. 
Consequently the scaling of the FA model in d = 3 is Gaussian. 

In Ref. |4|, the authors found that the directed percolation (DP) fixed point 
is relevant to the FA model in three dimensions and at low temperatures. This 
conclusion seems unsatisfactory: we find no evidence for fluctuation corrections in our 
simulation data. Further, and more importantly, the DP fixed point is characterised by 
a diverging static lengthscale ^: for a process in the DP universality class, one expects 
(nj(0)nj(0))DP — ('^)dp be a scaling function of \\x{i) — x{j)\\/C, where ^ diverges 
at criticality. On the other hand, detailed balance with respect to the non-interacting 
energy function ((1} tells us that {ni{0)nj{0))FA — (^)fa ^ij ^^e FA model. The 
absence of this diverging lengthscale in the FA model indicates that the underlying 
physics is different from that of the DP fixed point. Finally, we note that the exact 
result (3 = 1 arises naturally from our RG treatment: this situation would appear more 
satisactory than the argument of Ref. |1] that the exponents z and u should be given by 
their DP values while P is independently fixed to the non-DP value /3 = 1 by detailed 
balance. 

We argue that the conclusions of Ref. ^ regarding the DP fixed point are artefacts 
of an RG treatment that does not respect the presence of detailed balance and of the 
duality symmetry of the FA model; see also the comments after Equation ()78p. To 
remedy this, we have shown explicitly that writing the action as in ()79|1 allows one 
to perform the RG analysis in a way that preserves these symmetries. We attribute 
the apparent fluctuation effects in the data for the relaxation time of Ref. ji], which 
were derived from the persistence function P{t), to a combination of pre-asymptotic 
corrections in c and possible finite size effects. It should be emphasised that pre- 
asymptotic corrections are substantial: our simulations for X2{t) show that branching 
rates of c = 10~^ are still too large to see the true critical scaling. Similar comments 
apply to the persistence function: we do find Gaussian scaling as shown in Figure El but 
to see this clearly requires very small c. Note that in order to obtain data for c = 10~^ 
and c = 10~^ we used virtual simulation boxes of size 640^ and 2000^, respectively J. 
These cannot be reached by conventional lattice-based codes so that the data of Ref. |lj 
were of necessity taken from smaller systems, with potentially significant finite size 
effects. Although the DP exponents zbp = 1.9 and z/dp = 0.58 in li = 3 [H] differ 

I While our coordinate-based MC algorithm is extremely memory efficient (in a 2000^ lattice and at 
c = 10^^ there are on average only 8000 defects) keeping track of the persistence status is a problem 
for large lattices. We associate single bits with the persistence status of each lattice site so that we can 
track persistence in lattices up to 1000'^ using a memory block of 125Mb. At c = 10~^ we split the 
2000^^ virtual simulation box into eight 1000'^ blocks but only track persistence in four of them to limit 
the memory requirement to 500Mb. 
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Figure 5. Simulation data (symbols) for the equilibrium persistence function P{t) 
in the d — 2> FA model with branching rates c = 10~^, 10"'*, 10"^ and 10"^. The main 
panel shows data collapse for Gaussian scaling of the relaxation time r ^ 
inset demonstrates that the data are inconsistent with DP scaling tdp ^ c 
used the coordinate-based continuous- time MC algorithm to measure Pit) in direct 
equilibrium simulations; since P{t) is effectively a correlation function of events at a 
continuous range of times, the duality relations for two-time correlations cannot be 
used. The size of the virtual simulation box is 64'^, 200^, 640'^ and 2000^, respectively, 
which should be sufhcient to avoid finite-size effects. Results are averaged over 
1000, 100, 10 and 3 repeats, again in order of decreasing c. We expect relative errors 
of no more than 1% in the data shown. 



only slightly from the Gaussian ones, our data clearly allow us to rule them out: the 
dashed line in Figure El represents DP scaling of the dynamical correlation length and 
is inconsistent with our data. The inset in Figure El demonstrates similarly that the 
persistence functions P{t) do not collapse when plotted against t/rop with the DP 
scaling tdp ~ a rather similar picture - thus not shown - is obtained when 

plotting the two-point susceptibility X2(^), Figure El against t/rop. 

We conclude that the DP fixed point is irrelevant for the FA model, for the same 
reason that it is irrelevant for the AA model and for parity conserving models of 
branching and annihilating random walks: these models possess extra symmetries that 
must be preserved in RG calculations and, in the case of the FA model, lead to Gaussian 
scaling. This conclusion parallels that of Cardy and Tauber who showed that an 
early paper ^Hj on the parity conserving react ion- diffusion system {A SA, 2A 0) 
had produced a similar erroneous conclusion that the upper critical dimension was four 
and the exponents those of DP. 
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5. Conclusion 

To summarise, we showed in Section El that the FA and AA models with hard core 
exclusion share the same correlation functions (at equilibrium, and considering a single 
time-difference). This was established by means of an exact mapping at the level of 
the master equation. An important generalisation which includes additional diffusive 
processes demonstrated that the same mapping connects more generally the reaction- 
diffusion models with reversible coagulation A + A A and reversible annihilation 
A + A ^ 0. Further, we showed in Section El that the bosonic versions of the FA and 
AA models are appropriate effective theories for the low temperature limits of the hard 
core models and have analogous symmetries and relations between each other. Finally, in 
SectionElwe discussed the critical properties of the bosonic models using renormalisation 
group arguments. Implementing the mapping at the level of the field-theoretic action 
showed that the FA (and more generally A + A ^ A) model renormalises like the AA 
(or A + A 0) model. We find that the directed percolation fixed point is irrelevant to 
the FA model, because of the presence of detailed balance and of an additional hidden 
symmetry inherited from the parity symmetry of the AA model. Instead, the A+A^A 
model and its special case, the FA model, have upper critical dimension two; detailed 
balance together with the two exactly known scaling exponents is sufficient to find all 
exponents exactly also in < 2. 

From the point of view of reaction-diffusion systems and, more generally, non- 
equilibrium stochastic models, the most significant outcome of this work is the result 
that a hidden symmetry suppresses fluctuations in the A + A ^ A model and lowers its 
upper critical dimension to dc = 2. The mapping to A + A ^ can, however, also be 
used to more quantitative purposes. For example, it enables one to calculate new exact 
results for two-time non-equilibrium correlation and response functions in li = 1; we will 
report on these shortly jTHj. The results of such an analysis are instructive also more 
generally with regard to non-equilibrium fluctuation-dissipation relations for activated 
dynamics pi . In fact the FA model is an almost paradigmatic example of such dynamics, 
given that any evolution away from a metastable state containing only isolated defects 
requires the thermal excitation of additional defects. 

From a different angle, one may ask what our results have to say about the 
usefulness of the FA model for capturing the qualitative behaviour of structural 
glasses [3 El El- We have seen that in the physically relevant case of three spatial 
dimensions, fluctuation effects at low defect densities are of a classical (mean-fleld, 
Gaussian) nature. Nevertheless, the models will still exhibit a degree of dynamical 
heterogeneity. Violations of the Stokes-Einstein relation jS] may also persist, but 
will be at most by a constant (rather than diverging) factor as c — > 0; this is 
consistent with simulation results fHl- In summary, 'glassy' effects will be present, 
but probably rather weak. This is consistent with the fact that FA models also have 
relatively benign, Arrhenius-type increases of relaxation time scales at low temperature: 
T ~ ~ exp(2/T) in (i > 2 as we saw above. These models are therefore suitable 
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at best for modelling for what are known as strong glasses. For fragile glasses with 
their super-Arrhenius timescale divergences, models with facilitation by more than one 
spin - or with directed constraints - will inevitably have to be used. Their much 
more cooperative dynamics j2] continues to make them physically attractive models for 
understanding non-trivial aspects of glassy dynamics. 
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Appendix A. Large-S* expansion 



To see mathematically the equivalence between the bosonic and hard core models in the 
limit of small particle densities, one can replace the spin-1/2 operators of the hard core 
case by their spin-S* analogues and perform a formal large-S* expansion, valid for states 
with small density. For example, we can define a new Master operator by generalising 
£fa from (fTH|l to S* > 1/2, as follows: 
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{St - V2S) 



S Sj^ 
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{S, 



2Sc) + {i^ j) 



(A.l) 



where Sf etc are operators in the spin-S* algebra. In the spin-half case we have S = 1/2, 
{S + Sj)/{2S) = s^sj, {S — SI)/{2S) = s~'sf and so recover immediately £sfa = -^fa- 
The Master operator £sfa describes a system in which the number of particles on 
each site is restricted to the range < Ui < 2S; the particle number operators are 
hi = Sf + S. Configurations {ui} are mapped onto kets YYiiSf / \/2S)'"''\^) where |0) 
is the empty state as before; since nj|0) = is equivalent to |0) = —5*10) this state 
is, in spin language, fully polarised in the {—z) direction. Probabilities for transitions 
between states in some time interval t are then given by 



P{n[}^[n.}{t) = (0| 
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where we have introduced the coefficients F^^^ = {2SY"'{Q\{SfY{SfY\Q) for ease of 
writing. These obey the recursion F5^„+i = Fs„(n + 1)[1 — n/{2S)], yielding exphcitly 
r^.n = n\{2S)\ /[{2S — n)!(2S')"]. From ()A.2jl . conservation of probability requires that 
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which can be verified by direct calculation. We identify this left eigenstate as the 
projection state; it is analogous to (e| and (e| for hard core and bosonic models, 
respectively. 
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Using S+(5+/V2S)"|0) = V2S{S+/V2S)''+^\0) and 5r(S+/V2^)"|0) = V2S x 
n[l — n/ (2S)]{S^ / V2Sy''^^\0) one easily checks that the microscopic rates in the model 
defined by £sfa are 

riiUj {ui + l)nj, rate cnj\l - 72^/(25')], ,^ 

{ni + l)nj UiUj, rate nj{ni + — UiJ (23)]'^ 

These obey detailed balance with respect to the stationary state 

We have defined a family of interpolating models with increasing 5* that allow us to 
gradually remove the hard core constraints. We will now use a large-S* expansion to show 
that the models without constraints coincide with the bosonic models defined above. As 
long as there is no qualitative change in behaviour on increasing S we therefore expect 
the bosonic models to be suitable effective theories for the low temperature (small c) 
behaviour of the hard core ones. An example of a qualitative change that would render 
the large-S* expansion invalid is a transition to a quantum disordered state as 5* is 
reduced PU|: there is clearly no such singular behaviour here. Indeed, in our case the 
stationary states ()A.5|) of the models are known for general S. Bearing in mind that 
Tsfi = Ts,i = 1, they are all of effectively the same form if c is small so that only states 
with Ui = 0,1 have significant probability. 

Our claim that the above interpolating model becomes equivalent to the bosonic one 
in the limit S —>■ oo can be confirmed directly from ()A.4|1 : as long as c is small so that 
the relevant values of Ui stay small compared to 25*, the large-S* limit of the transition 
rates gives the bosonic model (jH). Correspondingly, the stationary state ()A.5jl becomes 
the bosonic one in this limit since — n!. 

More formally, one can establish the large-S* limit of our interpolating model by 
looking at the Hermitian version of the Liouvillian. Using detailed balance this is defined 
as HsFA = [n^c-"^/']/:sFA[niC"-/'] or exphcitly 



HsFA — 



(A.6) 



Then we can use the Holstein-Primakov representation 

S^ = ala^-S, S+ = al{2S ~ alaiY^^ (A.7) 

and take the large-S" limit by approximating 2S — ajai = 2S — fii ^ 2S everywhere. This 
assumes again that c is small enough so that all relevant states have particle numbers 
rii <^ 2S at each site. In spin language, all the states of interest are then localised on 
a small part of the surface of the spin sphere, and the non-trivial structure of the spin 
algebra can be neglected in favour of a simple bosonic one via ^ V2Sal Taking the 
S ^ oo limit as explained, we get 

^^SFA ^ J][4a,(aJ - Vc){ai - V^) + {i ^ j)]. (A.8) 



Mappings between reaction- diffusion and kinetically constrained systems 



32 



This coincides with the Hermitian form (j49|) of the Liouvilhan of the bosonic FA model 
as claimed, with the expected correspondence c = c. An exactly analogous procedure 
can be applied to construct a family of models that interpolates smoothly between the 
hard core and bosonic AA models. We therefore expect that the bosonic FA and AA 
models will be appropriate effective theories for their hard core counterparts at small 
particle densities. 

Appendix B. Random walk survival probabilities 

For the scaling analysis of equilibrium correlation functions given in Section 14.31 we 
required particular random walk survival probabilities; these are derived in the following. 
Consider a pair of diffusing defects. We can think of these as random walkers on a d- 
dimensional hypercubic lattice Z*^; whenever they occupy NN sites, where their positions 
have distance \\x2 — = 1, they coagulate instantaneously. We are interested in the 
probability Vd{x) that the walkers survive to infinite time, i.e. never coagulate. This 
survival probability depends on the spatial dimensionality d and on the initial separation 
X = X2 — Xi of the walkers. The distance vector X2 — Xi also performs a random walk, 
with twice the effective diffusion constant. The problem is therefore to calculate the 
probability Vd{x) that a random walker starting from position x will never reach one of 
the NN sites of the origin. If we picture these sites as absorbing, then Qd{x) = 1 —Vd{x) 
is the probability that the walker is absorbed eventually; for absorption at the origin 
itself, these quantities are well known. 

The key insight is that, in its first step, the walker randomly moves to one of the 
NN sites of x; we write these as y G N{x). The absorption probability starting from 
X is therefore the average of those that would be obtained when starting from any of 
these NN sites: 



The only exception to this relation is the case where x itself is already an absorption 
site so that Qd{x) = 1. We can correct for this by adding a source term at these sites; 
the latter then has to be set at the end of the calculation to give the correct values of 
Qd{x) at the absorption sites. Since all lattice directions are equivalent, all 2d source 
terms will be equal and we can write 




(B.l) 
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For the Fourier components Qk = Qd{x)e ^'"'^ this gives 




(B.3) 
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where a = 1, . . . ,d labels the lattice directions. The second, integral form of the result 
makes the reverse Fourier transform simple: 



Qdix) 



POO 

-6^,0 + d / dre-^^TT/,. 
Jo 



(B.4) 



where the /„ are modified Bessel functions, /n(T) = ^ /o'^dA; cos(nA;)e'^™'^^. Choosing 
V to ensure that Qd{x) = 1 for a; G A^(0) then gives for x 0, 

/„°°dre-'^%^-i(r)/i(r) 

We remark that for x = one finds, by retaining the Sxfi term, that Qd{0) = 1 as 
it must be: starting from the origin, already the first move reaches an absorption site. 
Given the derivation of ()B.5j) . it is not too surprising that the result is similar to that 
for the standard case of absorption at the origin, where one finds the same expression 
but with the denominator replaced by J dr e~'^'^/Q (r). 

In dimensions d = 1,2 the original pair of walkers always coagulate eventually, 
regardless of the initial separation and thus Vd{x) = 0. (This can be seen formally from 
the fact that both the numerator and denominator integrals in ()B.5|1 are dominated by 
their divergent tails; these have a;-independent prefactors, giving Qd{x) = 1.) In ci > 2, 
on the other hand. Equation ()B.5|1 yields a nonzero probability Vd{x) that the walkers 
survive indefinitely. The two particular values we need in the main text are 

rdre-3-/2(r)/2(r) , , 

Pa 2, 0,0=1- tW-^ 3.r2 , ^ 0-50166, B.6 

Jo dre ^-/o'(r)/i(r) 

^3(1, 1, 0) = 1 - ° ; ^ 0.35872. (B.7) 

Jo dre-^-/olr)/i(r) 

The rather significant difference between the numerical values fjB.6|) and ()B.7|) has a 
simple reason: when starting at initial separation x = (1, 1, 0) there are 12 possible first 
moves since each walker has 6 NN sites to move to. Out of these, 4 lead to the walkers 
being on NN sites where they coagulate instantaneously, while for x = (2, 0, 0) only 2 
moves produce this outcome. There are similar differences in subsequent moves which 
accumulate to the numbers given above. 
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